\documentclass[10pt,a4paper]{article} 

\usepackage{ctex} % 中文支持
\usepackage[top=2.5cm, bottom=2.5cm, left=2.5cm, right=2.5cm]{geometry} % 页边距
\usepackage{amsmath, amssymb} % 数学公式与符号
\usepackage{graphicx}
\usepackage{pythonhighlight}
\usepackage{url} 

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\usepackage{titling}
\setlength{\droptitle}{-2cm} % 标题上移

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%文档的题目、作者与日期
\author{五六七 }
\title{合金强度与回归分析 }

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
\begin{document}

\maketitle

\begin{abstract}
使用回归分析的方法来研究合金的强度与碳含量的关系。
\end{abstract}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\setcounter{tocdepth}{2}
%\renewcommand\contentsname{目录}
%
%\renewcommand {\baselinestretch} {1.3}\normalsize 
%\tableofcontents 
%\renewcommand {\baselinestretch} {1.0}\normalsize


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{问题描述}
合金的强度与其中的碳含量有密切的关系。现在从生产中收集到一些数据，如表所示，
其中 $x$ 为碳的含量，$y$ 为合金的强度。

\begin{table}[ht!]\centering
\caption{合金强度与碳含量的观测数据} \vspace{0.2cm}
\begin{tabular}{|c|c|c|c|c|c|c|c|c|c|c|c|c|}\hline
$x(\%)$ &0.10 &0.11 &0.12 &0.13 &0.14 &0.15 &0.16 &0.17 &0.18 &0.20 &0.22 &0.24 \\ \hline 
$y(kgf/mm^2)$ &42.0 &42.5 &45.0 &45.5 &45.0 &47.5 &49.0 &51.0 &50.0 &55.0 &57.5 &59.5 \\ \hline 
\end{tabular}
\end{table}

拟合一个函数关系式 $y=f(x)$, 并检验结果。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{回归分析基本概念}

函数拟合方法是根据样本的观测值，建立一个拟合数据的函数表达式，并估计拟合函数中的参数值的方法。
回归分析考虑的几个问题：
\begin{enumerate}
\item  建立因变量与自变量（一个或多个）之间的经验公式，也称为回归模型。
\item  根据样本数据，使用最小二乘法，估计经验公式中的参数。
\item  检验每个自变量在回归模型中是否显著。
\item  检验回归模型是否显著。
\item  使用回归模型对因变量进行预测。
\item  调整自变量对因变量进行控制。
\end{enumerate}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{最小二乘法的基础知识 }

设有观测数据 $(x_i,y_i), i=1,2,\cdots, n$. 最小二乘法是求参数 $\beta_0, \beta_1$ 使得
\begin{eqnarray}
Q = \sum\limits_{i=1}^n (y_i - \beta_0 - \beta_1 x_i)^2 
\end{eqnarray}
达到最小值。计算 $Q$ 关于 $\beta_0, \beta_1$ 的偏导数，可得
\begin{eqnarray}
\frac{\partial Q}{\partial \beta_0} &=& \sum\limits_{i=1}^n 2(y_i - \beta_0 - \beta_1 x_i)(-1), \\ 
\frac{\partial Q}{\partial \beta_1} &=&  \sum\limits_{i=1}^n 2(y_i - \beta_0 - \beta_1 x_i)(-x_i). 
\end{eqnarray}

当 $Q$ 达到最小值时，这两个偏导数都等于零。因此有
\begin{eqnarray}
n\beta_0 + (x_1+\cdots+x_n)\beta_1 &=& y_1+\cdots+y_n, \\ 
(x_1+\cdots+x_n)\beta_0 + (x_1^2+\cdots+x_n^2)\beta_1 &=& x_1y_1+\cdots+x_ny_n. 
\end{eqnarray}
由此可以求出 
\begin{eqnarray}
\begin{bmatrix}\hat{\beta}_0 \\ \hat{\beta}_1 \end{bmatrix} = (X^tX)^{-1}X^tY, 
\end{eqnarray}
其中 
\begin{eqnarray}
X= \begin{bmatrix} 1&x_1 \\ 1&x_2 \\ \vdots & \vdots \\ 1&x_n \\ \end{bmatrix}, \,\,\, 
Y= \begin{bmatrix} y_1 \\ y_2 \\ \vdots \\ y_n\end{bmatrix}. 
\end{eqnarray}


%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{建立模型}

设函数关系表达式为 
\begin{eqnarray}
y=\beta_0 + \beta_1 x + \varepsilon,
\end{eqnarray}
其中 $\beta_0, \beta_1$ 为待定参数，$\varepsilon$ 为随机误差。代入数据，可得12个方程，远多于参数个数。将这些方程写成矩阵形式，可得
\begin{eqnarray}
\begin{bmatrix}
42.0 \\ 42.5 \\ 45.0 \\ 45.5 \\ 45.0 \\ 47.5 \\ 49.0 \\ 51.0 \\ 50.0 \\ 55.0 \\ 57.5 \\ 59.5 \\ 
\end{bmatrix}
=
\begin{bmatrix}
1&0.10 \\ 1&0.11\\ 1&0.12\\ 1&0.13\\ 1&0.14\\ 1&0.15\\ 1&0.16\\ 1&0.17\\ 1&0.18\\ 1&0.20\\ 1&0.22\\ 1&0.24 \\ 
\end{bmatrix}
\cdot
\begin{bmatrix} \beta_0 \\ \beta_1 \end{bmatrix} 
+ 
\begin{bmatrix}
\varepsilon_1 \\ \varepsilon_2 \\ \varepsilon_3 \\ \varepsilon_4 \\ \varepsilon_5 \\ \varepsilon_6 \\ \varepsilon_7 \\ \varepsilon_8 \\ \varepsilon_9 \\ \varepsilon_{10} \\ \varepsilon_{11} \\ \varepsilon_{12} \\ 
\end{bmatrix}. 
\end{eqnarray}
简单记为
\begin{eqnarray}
Y=X\beta + \mathcal{E},
\label{eq-3}
\end{eqnarray}
其中 $Y$ 为应变量数据的列向量，$X$ 为常数系数1与自变量数据组成的列向量，$\beta$ 为参数列向量，$\mathcal{E}$ 为误差列向量。
根据最小二乘法的思路是，在等式 (\ref{eq-3}) 的两边同时左乘 $X^t$, 根据 $X^t\mathcal{E}\approx 0$, 可得 
\begin{eqnarray}
X^tY=X^tX\beta,
\end{eqnarray}
其中 $X^tX$ 是 $2\times 2$ 阶的可逆矩阵。因此可得
\begin{eqnarray}
\hat{\beta} = (X^tX)^{-1}X^tY. 
\end{eqnarray}
这就是最小二乘法的参数估计公式。

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{编程计算}

载入数值计算模块 numpy, 统计模型模块 statmodels, 和画图模块 matplotlib. 
\begin{python}
import numpy as np
import statsmodels.api as sm
import matplotlib.pyplot as plt
\end{python}

从文本文件载入数据，取第0行保存为自变量 x，取第1行保存为因变量 y. 
\begin{python}
d = np.loadtxt('data10_1.txt')
x0 = d[0]
y0 = d[1]
\end{python}

创建一个字典型数据作为最小二乘回归函数 ols 的数据参数，使用模型公式 y $\sim$ x. 拟合线性回归模型，结果保存为 re. 
\begin{python}
d = {'x':x0, 'y':y0}
re = sm.formula.ols('y~x', d).fit() 
\end{python}

输出回归模型的结果。
\begin{python}
print(re.summary())
\end{python}

模型输出结果如表所示。其中包含了很多信息，包括参数的估计值、每个参数的显著性检验、模型的显著性检验、残差的分布等。
\begin{python}
                            OLS Regression Results                            
==============================================================================
Dep. Variable:                      y   R-squared:                       0.977
Model:                            OLS   Adj. R-squared:                  0.975
Method:                 Least Squares   F-statistic:                     429.2
Date:                Sat, 21 Oct 2023   Prob (F-statistic):           1.52e-09
Time:                        21:54:50   Log-Likelihood:                -14.758
No. Observations:                  12   AIC:                             33.52
Df Residuals:                      10   BIC:                             34.49
Df Model:                           1                                         
Covariance Type:            nonrobust                                         
==============================================================================
                 coef    std err          t      P>|t|      [0.025      0.975]
------------------------------------------------------------------------------
Intercept     28.4835      1.030     27.648      0.000      26.188      30.779
x            129.0094      6.227     20.716      0.000     115.134     142.885
==============================================================================
Omnibus:                        3.342   Durbin-Watson:                   2.356
Prob(Omnibus):                  0.188   Jarque-Bera (JB):                1.843
Skew:                          -0.957   Prob(JB):                        0.398
Kurtosis:                       2.840   Cond. No.                         24.4
==============================================================================
\end{python}

输出已知数据的异常值的检验，输出残差的方差。
\begin{python}
print(re.outlier_test()) 
print(re.mse_resid)
\end{python}

输出的残差的方差如下。
\begin{python}
print(re.mse_resid)
0.8221698113207546
\end{python}

计算拟合值（模型的预测值），计算置信的下限和上限，计算置信半径。
\begin{python}
pre = re.get_prediction(d)
df = pre.summary_frame(alpha=0.05)
dfv = df.values 
low, upp = dfv[:,4:].T
r = (upp-low)/2 
num = np.arange(1, len(x0)+1)
\end{python}

画出散点图和拟合直线。
\begin{python}
fig=plt.Figure()
ax=fig.add_subplot(111)
ax.plot(x0,y0,'o')
ax.plot(x0,dfv[:,0],'r-')
ax.plot(x0,dfv[:,4],'b--')
ax.plot(x0,dfv[:,5],'b--')
\end{python}

\begin{figure}[ht!]\centering
\includegraphics [height=4cm, width=7cm]{alloy_regression_line.png}
\caption{散点图和回归直线 }
\end{figure}

画出残差图。
\begin{python}
plt.errorbar(num, re.resid, r, fmt='o')
\end{python}

\begin{figure}[ht!]\centering
\includegraphics [height=4cm, width=7cm]{alloy_error_bar.png}
\caption{残差图 }
\end{figure}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{检验模型}
从回归结果的表中可以看出，两个参数的p值约等于零，因此都是显著的。模型的$R^2$ 等于0.977, 因此模型也是显著的。



%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\section{回答问题}

根据一元线性回归分析，合金强度与含碳量呈现正向的线性关系，拟合函数的表达式为 
\begin{eqnarray}
\hat{y} = 28.4835 + 129.0094 x, 
\end{eqnarray}
这也称为经验回归方程。误差项的方差为 
\begin{eqnarray}
\hat{\sigma}^2 = 0.8222. 
\end{eqnarray}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%\section{参考文献 }
\begin{thebibliography}{99}

%\bibitem{dingtongren} 丁同仁、李承治，常微分方程教程，高等教育出版社，2022年3月第三版。
\bibitem{sishoukui-2} 司守奎,孙玺菁. \emph{Python数学建模算法与应用}, 国防工业出版社. 2022年1月第1版. 
\bibitem{hexiaoqun-ara} 何晓群. \emph{应用回归分析(R语言版)}. 电子工业出版社. 2017年7月第1版. 
\bibitem{maosisong} 茆诗松, 程依明, 濮晓龙. \emph{概率论与数理统计教程}. 高等教育出版社. 2019年11月第3版.


\end{thebibliography}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%

\end{document}

%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%
%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%%





